#include <thrust/gather.h>
#include <thrust/device_vector.h>
#include <thrust/scan.h>
#include <thrust/device_ptr.h>
#include "mian_util.h"
#include "result.h"
#include <cutil_inline.h>
struct resultCmp {
	__host__ __device__
		bool operator()(const h_result& o1, const h_result& o2) {
			return o1.pos < o2.pos;
		}
};

struct isZero {
	__host__ __device__
		bool operator()(const uint32 o1) {
			return o1 == 0;
		}
};

__device__ uint8 bit_table1[16] = {0,0,0,0,0,0,0,0,1,2,4,8,0x10,0x20,0x40,0x80}; 

__device__ void putBit(uint8* dest, uint32 pos, bool bit) {
	int curByte = pos >> 3;
	int curBit = pos & 7;
	dest[curByte] |= bit_table1[(uint32)bit<<3 + curBit];
}
/*
__device__ bool getBit(uint8* source, uint64 pos) {
	int curByte = pos/8;
	uint8 tempByte = source[curByte];
	tempByte >>= (7- (pos%8));
	bool bit = tempByte & 1;
	return bit;
}
*/
//inline __device__ uint32 curBytePos(x){ return (x)/8; } 
//data.x = dataCode
//data.y = dataLen

__device__ uint32 d_golomb_coding_int(uint32 num, uint32 q, uint32 r, uint8* dest, uint32 pos)
{
	uint32 j;
	uint32 len = 0;
	for(j=0; j<q; j++) {
		putBit(dest, pos + len, true);
		len++;
	}
	putBit(dest, pos + len, false);
	len++;

	int v = 1;
	for (j = 0 ; j < r; j++)
	{            
		putBit(dest, pos + len, v & num );
		len++;
		v = v << 1;         
	}
	return len; 
}
__device__ bool d_golomb_coding_4int(uint4* d_data, uint8* outBuffer, uint64 curPos, uint32 r)
{
	uint32 curBytePos = curPos >> 3;
	uint4 data = ((uint4*)d_data)[0];
	uint8* out = outBuffer + curBytePos;
	uint32 pos = curPos & 7;
	uint32 q1,q2,q3,q4;

	uint32 num1 = data.x;
	uint32 num2 = data.y;
	uint32 num3 = data.z;
	uint32 num4 = data.w;

	q1 = num1 >> r;
	q2 = num2 >> r;
	q3 = num3 >> r;
	q4 = num4 >> r;

	pos += d_golomb_coding_int(num1, q1, r, out, pos);

	pos += d_golomb_coding_int(num2, q2, r, out, pos);
	
	pos += d_golomb_coding_int(num3, q3, r, out, pos);

	pos += d_golomb_coding_int(num4, q4, r, out, pos);
	
	return true;
}
__device__ void d_golomb_coding_copy(uint8* out, uint8* output, uint64 lower, uint64 upper)
{
	uint32 byte_lower = lower >> 3;
	uint32 byte_upper = ceil((double)upper / 8.0);
	uint32 i = 0;
	for(i=byte_lower; i<byte_upper; i++)
		output[i] |= out[i - byte_lower];
	
}
//encoding the data set by even then odd ones to avoid write conflict
__global__ void d_golomb_coding_gpu_2pass(uint32* data, uint32 data_size, uint64* bitPos, uint8* output, uint32 r)
{
	uint32 i;
	uint4* data4 = (uint4 *)data;
	uint32 data_size4 = ceil((double)data_size / 4.0); 
	i = blockDim.x * blockIdx.x  + threadIdx.x;
	__shared__ uint8 out[128 * 4]; 
	uint32 pos_offset = i * 4;
	while((i * 2) < data_size4) {
		pos_offset = i * 4;		
		d_golomb_coding_4int(data4 + i*2, out, bitPos[pos_offset*2], r);
		pos_offset += 4;
		d_golomb_coding_4int(data4 + i*2+1, out, bitPos[pos_offset], r);
		d_golomb_coding_copy(out, output, bitPos[pos_offset - 4], bitPos[pos_offset + 4]);
		i += blockDim.x * blockIdx.x;
	}
}
/*
__device__ bool d_golomb_coding_int_buffer(uint4* data, uint8* outBuffer, uint32 bufferLenPerInt, uint4* dataLen, uint32 M)
{
	uint8 out[bufferLenPerInt * 4];
	int r = log((double)M)/log((double)2);
	uint4 data = ((uint4*)data)[0];
	uint4 len;
	uint32 pos = 0;
	uint32 j,q1,q2,q3,q4;

	uint32 num1 = data.x;
	uint32 num2 = data.y;
	uint32 num3 = data.z;
	uint32 num4 = data.w;

	q1 = num1 / M;
	q2 = num2 / M;
	q3 = num3 / M;
	q4 = num4 / M;

	len.x = q1 + 1 + r;
	len.y = q2 + 1 + r;
	len.z = q3 + 1 + r;
	len.w = q4 + 1 + r;

	if(len.x + len.y + len.z + len.w > 4 * bufferLenPerInt)
		return false;
	d_golomb_coding_int(num1, q1, r, out, pos);
	pos += len.x;

	d_golomb_coding_int(num2, q2, r, out, pos);
	pos += len.y;
	
	d_golomb_coding_int(num3, q3, r, out, pos);
	pos += len.z;

	d_golomb_coding_int(num4, q4, r, out, pos);
	pos += len.w;
	
	for(j=0; j<bufferLenPerInt; j++)
		((uchar4*)(outBuffer))[j] = ((uchar4*)out)[j];
	dataLen[0] = len;
}

__device__ bool d_golomb_coding_thread(uint4* data, uint8* buffer, uint4* dataLen, uint64* curPos, 
		uint32 maxBufferLen, uint32 M)
{

	uint32 j,q1,q2,q3,q4;
	uint64 pos = curPos;
	int r = log((double)M)/log((double)2);
	uint4 data = ((uint4*)data)[0];
	uint4 len;
	uint32 num1 = data.x;
	uint32 num2 = data.y;
	uint32 num3 = data.z;
	uint32 num4 = data.w;

	q1 = num1 / M;
	q2 = num2 / M;
	q3 = num3 / M;
	q4 = num4 / M;

	if(curBytePos( (q1 + 1) + (q2 + 1) + (q3 + 1) + (q4 + 1) + (r * 4) + pos) > maxBufferLen) {
		(*curPos) = 0;
		return false;
	}
	
	d_golomb_coding_int(num1, q1, r, buffer, pos);
	pos += q1 + 1 + r;
	len.x = q1 + 1 + r;

	d_golomb_coding_int(num2, q2, r, buffer, pos);
	pos += q2 + 1 + r;
	len.y = q2 + 1 + r;

	d_golomb_coding_int(num3, q3, r, buffer, pos);
	pos += q3 + 1 + r;
	len.z = q3 + 1 + r;

	d_golomb_coding_int(num4, q4, r, buffer, pos);
	pos += q4 + 1 + r;
	len.w = q4 + 1 + r;

	(*curPos) = pos;
	dataLen[0] = len;
}


__global__ void d_golomb_coding(uint32* dataBuffer, uint32 dataSize, uint32* dataLenBuffer,
		uint8* outDataBuffer, uint32 outDataBufferSize, uint32 M) {

	uint32 i = blockIdx.x * blockDim.x + threadIdx.x;
	uint32 t_id = threadIdx.x;
	uint4* data = (uint4*)dataBuffer;
	uint4* data_len = (uint4*) dataLenBuffer;
	const int max_buffer_size = outDataBufferSize / (blockDim.x * gridDim.x);
	const int out_buffer_size_per_int = outDataBufferSize / (dataSize * sizeof(int));
	const int out_buffer_size_per_int4 = out_buffer_size_per_int * 4;
	uint64 bitPos = 0;

	while(i  < dataSize){
//		d_golomb_coding_thread(data + i, outDataBuffer + t_id * thread_buffer_size, 
//					data_len + i, &bitPos, max_buffer_size, M);
		d_golomb_coding_int_buffer(data + i, outDataBuffer + out_buffer_size_per_int4*i, 
						out_buffer_size_per_int, data_len + i, M);
		i += blockIdx.x * blockDim.x;
	}		
}

void Glomb_collect_over_size_gpu(uint32* dataBuffer, uint32* dataLenBuffer, uint32* map, 
				uint32 dataSize, uint32* outDataBuffer, uint32* outSize)
{
	thrust::device_ptr t_values = device_ptr(dataBuffer);
	thrust::device_ptr t_map = device_ptr(map);
	thrust::device_ptr t_stencil = device_ptr(dataLenBuffer);
	thrust::device_ptr t_output = device_ptr(outDataBuffer);

	thrust::gather_if(t_map, t_map + dataSize,
			t_stencil,
			t_values,
			t_output,
			is_zero());
	(*outSize) = thrust::count(t_stencil, t_stencil+dataSize, 0);
}
*/
/*
void Golomb_coding_gpu(uint32* d_data, uint32 d_data_size, uint8* output)
{
	uint32* d_overSizeDataBuffer, *map;
	uint8* d_data_small_output, d_data_large_output;
	uint32* M, *d_overSizeDataNum;
	uint32 overSizeDataNum;
	uint32 buffer_size_per_thread;

	cutilSafeCall(cudaMalloc(&d_overSizeDataBuffer, sizeof(int), d_dataSize));
	cutilSafeCall(cudaMemset(d_overSizeDataBuffer, 0, d_dataSize * sizeof(int)));
	cutilSafeCall(cudaMalloc(&d_data_small_output, 1, d_dataSize * sizeof(int)));
	cutilSafeCall(cudaMemset(d_data_small_output, 0, d_dataSize * sizeof(int)));
	cutilSafeCall(cudaMalloc(&d_data_large_output, 1, d_dataSize * sizeof(int)));
	cutilSafeCall(cudaMemset(d_data_large_output, 0, d_dataSize * sizeof(int)));
	cutilSafeCall(cudaMalloc(&M, sizeof(int), 1));
	cutilSafeCall(cudaMalloc(&d_overSizeDataNum, sizeof(int), 1));
	
	Golomb_choseM(d_data, M);
	d_golomb_coding(d_data, d_data_size, sizeof(int), d_dataLen, d_data_small_output, M);
	
	Golomb_collect_over_size_gpu(d_data, d_dataLen, map, 
			d_data_size, d_overSizeDataBuffer, d_overSizeDataNum);
	cutilSafeCall(cudaMemcpy(&overSizeDataNum, d_overSizeDataNum, sizeof(int), 1, cudaMemcpyDeviceToHost));

	uint32 newBufferSize = d_dataSize * sizeof(int) / overSizeDataNum;
	Golomb_Coding_gpu(d_data, d_data_size, newBufferSize ,d_dataLen, d_data_large_output, M);

	Golomb_Merge_gpu(d_data_small_output, d_data_large_output, d_dataLen);
}	
*/
uint32 get_golomb_para_gpu(uint32* data, uint32 dataSize)
{
	thrust::plus<uint32> binary_op;
	thrust::device_ptr<uint32> t_data(data);
	uint32 init = 0;
	uint64 sum = thrust::reduce(t_data, t_data + dataSize, init, binary_op);
	
	return  (uint32)(sum / (uint64)dataSize);	
}

__device__ void d_golomb_collect_thread(uint8* d_output_buffer_ptr, uint64 lower, uint64 upper, uint8* d_out_ptr)
{
	uint32 lower_byte = lower / 8;
	uint32 upper_byte = upper / 8;
	if(lower_byte % 8 != 0)
		lower_byte += 1;
	if(upper_byte % 8 != 0)
		upper_byte -= 1;
	
	uint32 i;
	for(i=lower_byte; i<upper_byte; i++)
		d_out_ptr[i] = d_output_buffer_ptr[i - lower_byte];	
}

/*
uint8 all_one[9] = {0, 1, 3, 7, 0xF, 0x1F, 0x3F, 0x7F, 0xFF};

__device__ void d_golomb_collect_thread2(uint8* d_output_buffer_ptr1, uint64 d_output_buffer_bitlen1, uint8* d_output_buffer_ptr2, uint64 curPos, uint8* d_out_ptr)
{
	uint32 left = 8 - (curPos % 8);
	uint8 left_byte = d_output_buffer_ptr1[ (d_output_buffer_bitlen1 - left)/8];
	left_byte &= all_one[left] ;
	left_byte <<= 8 - left;
	d_out_ptr[0] = left_byte;

	uint32 right = curPos % 8;
	uint8 right_byte = d_output_buffer_ptr2[0];
	right_byte >>= left;
	d_out_ptr[0] |= right_byte;
	
}

__device__ void d_golomb_collect(uint8* d_output_buffer_ptr, uint32 d_buffer_size_per_int4, uint32 d_data_num,
				uint64* d_data_len_ptr, uint8* d_output_ptr)
{
	uint32 i = blockDim.x * blockIdx.x + threadIdx.x;
	uint64 curPos;
	uint32 d_buffer_per_thread = d_buffer_size_per_int4 * sizeof(int);
	while(i < d_data_num){
		d_golomb_collect_thread(d_output_buffer_ptr i * d_buffer_size_per_thread, 
					d_data_len_ptr + i * 4, d_data_len_ptr + (i+1)*4, d_out_ptr);
		i += blockDim.x * blockIdx.x;
	}
	i = blockDim.x * blockIdx.x + threadIdx.x + 1;
	uint64 lower, upper;
	while(i < (d_data_num/4)) {
		lower = (i - 1) * 4
		d_golomb_collect_thread2(d_output_buffer_ptr + (i-1)*4 );
	}
}
*/
__global__ void d_golomb_get_bit_pos(uint32* data, uint32 data_size, uint64* bitPos, uint32 r)
{
	uint32 i = blockDim.x * blockIdx.x + threadIdx.x;
	while(i < data_size) {
		bitPos[i] = (data[i] >> r) + 1 + r;
		i += blockDim.x * blockIdx.x;
	}
}

void golomb_coding_gpu_1pass(uint32* inBuffer, uint32 inBufferSize, uint8* outBuffer, uint32* outBufferSize)
{
	const int thread_num = 128, block_num = 512;
	uint32* d_data_ptr, *d_data_len_ptr;
	uint8* d_output_ptr, *d_output_buffer_ptr;
	
	const uint32 output_buffer_size_per_int = 4;
	const uint32 output_buffer_size_per_int4 = 4 * 4;
	uint32 total_output_buffer_size = inBufferSize * sizeof(int) * output_buffer_size_per_int;
	uint32 output_buffer_size_per_thread = total_output_buffer_size / (thread_num * block_num);

	cutilSafeCall(cudaMalloc(&d_data_ptr, sizeof(int) * inBufferSize));
	cutilSafeCall(cudaMemcpy(d_data_ptr, inBuffer, sizeof(int)*inBufferSize, cudaMemcpyHostToDevice));
/*	
	cutilSafeCall(cudaMalloc(&d_data_len_ptr, sizeof(int) * inBufferSize));
	cutilSafeCall(cudaMemset(d_data_len_ptr, 0, sizeof(int) * inBufferSize));
*/	
	cutilSafeCall(cudaMalloc(&d_output_ptr, sizeof(int) * inBufferSize));
	cutilSafeCall(cudaMemset(d_output_ptr, 0, sizeof(int) * inBufferSize));
	
//	cutilSafeCall(cudaMalloc(&d_output_buffer_ptr, total_output_buffer_size));
//	cutilSafeCall(cudaMemset(d_output_buffer_ptr, 0, total_output_buffer_size));

	uint32 M = get_golomb_para_gpu(d_data_ptr, inBufferSize);
	uint32 r = ceil(log((double)M) / log(2.0));
	M = pow(2.0, (double)r);
/*	
	d_golomb_coding<<<block_num, thread_num>>>(d_data_ptr, inBufferSize, d_data_len_ptr, d_output_buffer_ptr, 
						total_output_buffer_size, M);

	thrust::device_ptr t_data_len_ptr = device_ptr(d_data_len_ptr);
	thrust::exclusive_scan(t_data_len_ptr, t_data_len_ptr + inBufferSize);	
	d_golomb_collect<<<block_num, thread_num>>>(d_output_buffer_ptr, d_data_len_ptr, output_buffer_size_per_int4, d_output_ptr);
*/
	uint64* d_bitPos_ptr;
	cutilSafeCall(cudaMalloc(&d_bitPos_ptr, sizeof(uint64) * inBufferSize));
	cutilSafeCall(cudaMemset(d_bitPos_ptr, 0, sizeof(uint64) * inBufferSize));
	
	d_golomb_get_bit_pos<<<block_num, thread_num>>>(d_data_ptr, inBufferSize, d_bitPos_ptr, r);
	thrust::device_ptr<uint64> t_bitPos_ptr(d_bitPos_ptr);
	thrust::exclusive_scan(t_bitPos_ptr, t_bitPos_ptr + inBufferSize, t_bitPos_ptr, 0);	
	d_golomb_coding_gpu_2pass<<<block_num, thread_num>>>(d_data_ptr, inBufferSize, d_bitPos_ptr, d_output_ptr , r);
	
	uint64 out_size;
	cutilSafeCall(cudaMemcpy(outBuffer, d_output_ptr, sizeof(int) * inBufferSize,cudaMemcpyDeviceToHost));
	cutilSafeCall(cudaMemcpy(&out_size, d_bitPos_ptr+inBufferSize, sizeof(int), cudaMemcpyDeviceToHost));

	printf("after compress: %u\n", ceil(out_size / 8.0));

	cutilSafeCall(cudaFree(d_data_ptr));
	cutilSafeCall(cudaFree(d_output_ptr));
	cutilSafeCall(cudaFree(d_bitPos_ptr));
}
/*
int main(int argc, char* argv[])
{
	char* posFile = "posFile";	
	FILE* fin = fopen(posFile, "r");
	uint32 dataSize = 1000000;
	uint32* posData = (uint32*)malloc(dataSize * sizeof(int));
	uint8* posComp = (uint8*)malloc(dataSize * sizeof(int));
	uint32 compSize = 0;
	fread(posData, sizeof(int), dataSize, fin);
	golomb_coding_gpu_1pass(posData, dataSize, posComp, &compSize);

	fclose(fin);
	free(posData);
	free(posComp);
}
*/
